Effects of charge-dependent vibrational frequencies and anharmonicities in transport through 
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As a step towards a more realistic modeling of vibrations in single-molecule devices, we investigate the ef- 
fects of charge-dependent vibrational frequencies and anharmonic potentials on electronic transport. For weak 
phonon relaxation, we find that in both cases vibrational steps split into a multitude of substeps. This effectively 
leads to a bias-dependent broadening of vibrational features in current-voltage and conductance characteristics, 
which provides a fingerprint of nonequilibrium vibrations whenever other broadening mechanisms are sec- 
ondary. In the case of an asymmetric molecule-lead coupling, we observe that frequency differences can also 
cause negative differential conductance. 

PACS numbers: 73.23.Hk,73.63.-b,85.65.+h 
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Introduction. — Recently, several experiments have probed 
the electron-phonon coupling in single-molecule devices by 
detecting vibrational features in current-voltage characteris- 
tics (IVs)}'^'^''^ This coupling and its consequences for elec- 
tronic transport through molecules have also been at the focus 
of theoretical studies of /y£^ '^ii?'f'?4" and noise spectra,'"'" 
and it is perceived as a possible avenue towards the design of 
molecular devices. 

The signatures of electron-phonon coupling differ depend- 
ing on the transport regimes. Vibrational steps appear in IVs 
in the sequential-tunneling regime,^-^ "' " and in the conduc- 
tance dl/dV in the inelastic-cotunneling regime.** The visi- 
bility of these steps generally depends on the step heights and 
spacings, as well as their broadening. Three main mechanisms 
for broadening of vibrational steps in IV and conductance 
of single-molecule devices have been identified: Broadening 
induced by (i) temperature, (ii) vibrational dissipation,^ and 
(iii) molecule-lead tunneling.^ In this paper, we discuss an ad- 
ditional broadening mechanism which arises when going be- 
yond the one-mode harmonic approximation of previous mod- 
els. 

Essentially, all studies to date, e.g. Ref s . '5"6'7"8" 1 Ol 111 are 
based on two simplifying approximations. First, they restrict 
the description of the molecular vibrations to the harmonic ap- 
proximation, and second, they assume the vibrational frequen- 
cies to be identical for the different molecular charge states 
relevant for the transport. In devices with symmetric voltage 
splitting, the combination of these two approximations leads 
to strictly equidistant steps in IV for specific gate voltages. 

Here, we take a step towards a more realistic modeling of 
the vibrations by investigating both the case of vibrational 
frequencies depending on the molecular charge state, and 
the case of anharmonic oscillations within a Morse-potential 
model. Using a rate-equations approach valid for the regime 
of weak molecule-lead coupling, we calculate IVs and con- 
ductances, and show that strictly equidistant vibrational steps 
are indeed an artefact of the simplifications in previous mod- 
els. Remarkably, we find that the strength of direct vibrational 
relaxation can significantly alter the current-voltage charac- 
teristics. For strong relaxation, the extended models mainly 
lead to shifts of the step positions. In contrast, for weak vi- 
brational relaxation they effectively result in a bias-dependent 



broadening of vibrational steps due to the steps splitting into 
a multitude of closely spaced substeps. 

Model. — Our starting point is a generic model for a 
molecule coupled to metallic leads. Transport is 
assumed to be dominated by one spin-degenerate elec- 
tronic level with energy e (measured with respect to 
the zero-bias Fermi energies of the leads), tunable by a 
gate electrode. The system is described by the Hamil- 
tonian H = H^oi + -ff leads + ^^mix, where iJieads = 
Y.a=L,RY.p,cr'^pcipaCapa descrfbcs the non-intcracting 
leads, i/mix = Ea=L,i?,Ep,a {tacip^da +h.c.) the tunnel- 
ing between leads and molecule, and 



U P'^ 
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models the molecular degrees of freedom. Here, U is the 
charging energy for double occupation, (rft) annihilates 
(creates) an electron with spin projection cr on the molecule, 
and Tid = J2cr d]yda denotes the corresponding occupation- 
number operator Similarly, Capa ic\p^) annihilates (creates) 
an electron in lead a {a — L,R) with momentum p and spin 
projection cr. 

The kinetic and potential energy terms in Eq. Q re- 
fer to the nuclear motion. For definiteness, we consider a 
model with one dominant mode of vibrations. Then, X, 
P, and n denote the corresponding normal coordinate, mo- 
mentum, and reduced mass, respectively. For the singly 
and doubly charged molecular ions, the potential energy 
curves generally deviate from their neutral counterpart. In 
the spirit of the Born-Oppenheimer approximation, we take 
this into account by writing the potential energy in the form 
Vnj^{X) = J2n'''n{X) \n) {n\, using projectors onto the 
electronic ground states with fixed molecular charge n — 
0,1,2. 

The coupling between molecule and leads, parameterized 
by the tunneling matrix elements and t^j, is assumed to be 
weak, i.e. the tunneling-induced energy broadening 7 of elec- 
tronic levels is small, 7 ^ ksT, hojQ, where T denotes the 
temperature, and ojq the vibrational frequency. Then, i?,nix 
can be treated perturbatively and the solution of rate equa- 
tions is sufficient for current and noise calculations i i »i i -ii i " For 
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simplicity, we assume symmetric voltage splitting throughout 
this paper, and focus on the case of strong Coulomb blockade 
(U -> oo). 

In the absence of tunneling, the molecular states can be 
written as \n,q) with n denoting the charge state of the 
molecule and q the number of excited phonons. The rates 



L,R ^^q- 



for transitions 



q^q' ^a=L,H • ■ q^q' ■,a \n,q) — > 

I n', q' ) are calculated via Fermi's golden rule. The occupa- 
tion probabilities are determined by the rate equations 
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where the last term describes direct vibrational relaxation 
towards the equilibrium phonon distribution Pq"' — [1 
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^Q-Pqfi'^o ^jjjj j.^jg i^^^ ^Yhe steady-state current is 
calculated by solving the rate equations (|3 in the stationary 
case and evaluating 
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Recent papers investigating the influence of molecular vi- 
brations on current and noise have focused on the model 
Vn{X) — i/iWQ(X+\/2"-A£osc)^, i-S- vibrations are taken into 
account within the harmonic approximation and their frequen- 
cies are assumed to be independent of the molecular charge 
state. ^i^iiSiiSiii The position of the potential minimum explic- 
itly depends on the charge state, and its shift, measured in 
units of the oscillator length £o,sc = (^//^^^o)^^^. is character- 
ized by the electron-phonon coupling strength A. 

We extend this model as follows. In the first case ( IMU . we 
include frequency variations for the different charge states by 
using 



(Ml) 



In the second case iM2\ . we investigate the effects of anhar- 



monicities within a model based on the Morse potentially 



(M2) 



where D denotes the dissociation energy, (3 determines the 
inverse range of the potential, and ^osc = (V/^We)^/^ char- 
acterizes the range of the ground-state wave function with 
LUe = Pi'^D/ nY^^ . The number of bound states is given by 
+ 1 where j = ^2^D/I3h - 1/2. The Morse potential 
encompasses both bound and continuum states and is there- 
fore also suitable to describe dissociation processes.'^ In this 
paper, we restrict the discussion to situations in which transi- 
tions into continuum states can be neglected. 

An important qualitative difference compared to the har- 
monic approximation arises from the asymmetry of the Morse 
potential under parity transformations X — > ~X, which 
causes the overlap of vibrational wavefunctions and hence the 
transition rates to depend on the direction of the shift of the 
potential minima. In contrast to the case of symmetric po- 
tentials, the FC matrix elements now behave differently de- 
pending on the sign of the electron-phonon coupling A.'^ (For 
symmetric potentials the sign is irrelevant and A can be cho- 
sen to be a positive number) In the general case, this leads to 
a dependence of the stationary current on the sign of A. 

Calculations. — The calculation of the transition rates 
^q-^q' essentially reduces to the evaluation of Franck- 
Condon (FC) matrix elements,'** given by the overlap of the 

of the initial and final states 



oscillator wave functions 
\n,q) and 



^n.q 



dx (f>* Jx)(t)n'.,q'{x). 
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In harmonic models with loq = loi = 0J2 the FC matrix ele- 
ments are independent of n, n' and can be compactly written 
in terms of Laguerre polynomials, see e.g. Ref. J^, For the 
harmonic model ( IM1> with different vibrational frequencies 
ujn, and the Morse potential model ( IM2t . no such simple ex- 
pression is available^ and instead we evaluate the integrals 
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numerically, where Eqs. Q and (|6} correspond to (IMU and Since the transition rates WJ^^^ are bounded from above, 

( IM2> . respectively. Hq (x) denotes the Hermite polynomial of and the stationary occupation probabilities obey P^ for 

order q, L"{x) the generalized Laguerre polynomial, an — q ^ 00, the Eqs. (|2} and (|3} can be effectively solved for a fi- 

y^uj^/ujQ and a — e^5'3("'-")\/2A£„sc nite number of re levant phonon excitations. For the Morse po- 
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FIG. 1: (Color online) IV and conductance for a symmetric de- 
vice (t_L ~ in, symmetric voltage splitting) with slightly differ- 
ent vibrational frequencies for the neutral and charged molecule 
(wi = l.ltJo), and A = 1.2, U ^ oo. (a) /V at e = and 
ksT <^ huoo for strong and weak vibrational relaxation. While for 
equilibrated phonons, step locations are specified by Eq. 0, for un- 
equilibrated phonons a multitude of steps is generated, leading to an 
effective step-broadening even at low temperatures. For comparison, 
the inset shows results with identical frequencies uja = uoi. (b),(c) 
Conductance plots at ksT = 0.02fttJo showing the subsplitting and 
washing out of phonon features for unequilibrated phonons. 



tential, this treatment requires the restriction to a sufficiently 
high number of bound states and sufficiently low voltages so 
that all relevant excitations remain far below the dissociation 
threshold, q ^ [j] . 

Results and interpretation. — ^We first discuss results for 
model ( IM1> with harmonic vibrations and charge-dependent 
frequencies. Representative results are depicted in Fig. [2 For 
equilibrated phonons (r 0) and frequencies tjQ 7^ steps 
appear at different positions than in the identical-frequencies 
case, see Fig. 1(a). In principle, steps remain well-defined as 
shown in the corresponding conductance plot, see Fig. ^h), 
where only small shifts of the series of vibrational sidebands 
are caused by the entering of two different energy scales loq 
and LOi. In contrast, for unequilibrated phonons (r oo) a 
splitting of vibrational steps into a multitude of substeps is ob- 
served to cause an effective broadening, which increases with 
bias voltage, see Fig.[na),(c). We now give a more quantita- 
tive explanation for this behavior. 

For equilibrated phonons and temperatures 7 ^ fc^T <C 
fvjjQ^i, the phonon distribution Pq = P^' is strongly 
peaked for the phonon ground state q = Q. Consequently, 
phonon transitions predominantly occur in the channel ^ g. 
Energy conservation limits the possible excitation processes, 
and in the (eV, e)-plane the boundaries at which the transi- 
tions I 0, ) III?) and I 1, ) I 0, 9 ) become possible 
are given by 

e = Tey/2 + - El , and e - TeV/2 + E° - , (7) 

where E"^ = hwniq + 1/2) denotes the phonon energy of the 
state \n,q), and the lower (upper) sign refers to the left (right) 
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FIG. 2: (Color online) Current-voltage characteristics for an asym- 
metric device {t\ — 0.3ti, symmetric voltage splitting) in the un- 
equilibrated regime with \ = 1.1, ksT huo,U ~> 00. The com- 
bination of asymmetric coupling and differing vibrational frequen- 
cies is observed to lead to weak peaklike structures at the onsets of 
phonon steps with negative differential conductance regions (marked 
by arrows). 

lead. In addition, cuiTent only flows if electrons can traverse 
the molecule, leading to the condition 

\eV/2\>\e^E^ + E'o\, (8) 

which together with Eq. completely specifies the positions 
of all vibrational sidebands for equilibrated phonons. In par- 
ticular, even small frequency differences result in IVs with 
non-equidistant phonon step spacings,^^ and the central cross- 
ing point is now located at (eV, e) — (0, h[ijjQ — wi]/2) due to 
the deviation of the zero-point energies of the oscillators. 

In contrast, for unequilibrated phonons a multitude of ex- 
citation and de-excitation processes \ Q,q) and 
1 1, (7 ) I 0, g' ) are permitted, and Eq. must be replaced 
by 

e = TeV/2+ E°-El,, and e = TeV/2+E°-El. (9) 

However, the range of the contributing q, q' in Eq. (|9} now 
inherently depends on the nonequilibrium distribution P^, 
which makes concise statements about the occurence of cer- 
tain steps more difficult. Generally, larger bias voltages \eV\ 
cause the nonequilibrium distribution to widen, and thus in- 
crease the range of relevant q, q' in Eq. (|5Jl. This corresponds 
to a splitting of a vibrational step into an increasing number of 
closely spaced substeps as observed in Fig.^a),(c), resulting 
in an effective broadening of steps. 

Interestingly, for asymmetric devices (<l 7^ t r) we find that 
frequency differences between molecular charge states can in- 
duce negative differential conductance (NDC). A representa- 
tive example for this behavior is shown in Fig.|2] The split- 
ting of phonon steps causes peaklike structures at the onset 
of phonon steps. For certain parameters, such NDC features 
have been predicted for asymmetric coupling and harmonic 
vibrations with ~ uJi in Refs. i5|. Here, we report that 
qualitatively similar NDC behavior can arise from differing 
vibrational frequencies, for parameters where our numerical 
results do not exhibit NDC for luq = It is interesting to 
note that NDC at onsets of phonon steps is also observed in 
the experiment by Park ef oL ' '^ 

Results for the model ( IM2> with anharmonic vibrations are 
depicted in Fig.|3] While the main effects are similar to those 
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FIG. 3: (Color online) IV and conductance plots for a symmetric 
device (tL = tR, symmetric voltage splitting) with a Morse potential 
(j = 30), for A = 1.5, U ^ oo. (a) 11/ at e = and ksT < hcue 
for strong and weak vibrational relaxation. For equilibrated phonons 
steps along e = are located at voltages specified by Eq. 0, for 
unequilibrated phonons the step-splitting leads to an effective step- 
broadening even at low temperatures. (b),(c) Conductance plots at 
ksT — 0.02hLJe showing the subsplitting and broadening of phonon 
features for unequilibrated phonons. 

of model JMll l - steps in IV being non-equidistant, weak vi- 
brational relaxation leading to the splitting of steps into sub- 
steps - the underlying mechanism is different. In the Morse 
model ( IM2t . the potential curves are shifted depending on the 
molecular charge state, but are identical otherwise. Here, the 
appearance of different step spacings and splitting of steps is 
purely due to the fact that the energy spectrum of an anhar- 
monic oscillator is not equidistant. 

In analogy to model iMl\ . for equilibrated phonons the po- 
sitions of vibrational sidebands are fixed by substituting the 
eigenenergies of the Morse potential Eg = —D + huje{q + 
1/2) - huieXiq + 1/2)^ into Eqs. Q and (jS). Here, we have 
used X = (2i + 1)^- This leads to a conductance plot 
consisting of two sets of parallel lines with decreasing line 
distances, reflecting the decreasing spacing of eigenenergies 
of the Morse potential as the quantum number q increases, 
see Fig.|3la),(b). Since the phonon ground state energies are 
identical, Eq = Eq, the central crossing point is located at 



(ey,£) = (0,0). 

In the unequilibrated case, the anharmonicity of the poten- 
tial and the occurence of various phonon excitation and de- 
excitation processes give rise to a splitting of steps. As be- 
fore, the range of relevant q, q' in Eq. (|9} depends on bias 
and gate voltage, leading to a growing number of substeps 
with increasing voltage, and hence an effective washing out 
of phonon steps, see Fig.|3a),(c). 

Conclusions. — As a step towards a more realistic modeling 
of molecular vibrations and their consequences for electronic 
transport in single-molecule devices, we have investigated the 
effects of charge dependences of vibrational frequencies and 
anharmonicities in the sequential-tunneling regime. Even for 
small frequency differences and anharmonicities, we find that 
vibrational step spacings cease to be equidistant. In combi- 
nation with weak vibrational relaxation, both frequency dif- 
ferences and anharmonicities are shown to lead to a bias- 
dependent splitting of levels, effectively resulting in a broad- 
ening of phonon steps. For asymmetric molecule-lead cou- 
pling, we find that this mechanism can also lead to negative 
differential conductance behavior. We conclude that spac- 
ings of vibrational features in dl/dV provide information 
about charge dependence of vibrational frequencies and an- 
harmonicity of the potential. Whenever other broadening 
mechanisms play a secondary role, the bias-dependent sub- 
splitting or broadening acts as a fingerprint of nonequilibrium 
vibrations. 

It is interesting to note that the appearance of clusters of 
vibrational substeps due to the emergence of differing energy 
scales of phonon excitations resembles the splitting of tun- 
neling peaks observed in transport through ultrasmall metallic 
grains. In that context, Agam et al. argued that the Coulomb 
interaction effectively results in a bias-dependent splitting of 
resonance peaks in the nonequilibrium case, i.e. for suffi- 
ciently slow electronic relaxation. Finally, we remark that 
similar arguments also suggest a splitting of the conductance 
fine structure for magnetic single-molecule devices,^" when- 
ever the exchange coupling varies with the molecular charge 
state. 
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We point out that an analytical solution of the integral in Eq. J6) in 

terms of an alternating sum is possible.'^ However, for potentials 



with a large number of bound states the sum is intractable, and 
direct numerical quadrature is preferrable. 

We remark that model <M1> can still yield IVs with equidis- 
tant phonon steps for equilibrated phonons, if asymmetric voltage 
splitting (not considered here) is introduced. 



